Estimation of System Parameters in Discrete Dynamical Systems from Time Series 



O 

o 

(N 



p. Palaniyandi and M. Lakshmanaifl 

Centre for Nonlinear Dynamics, Department of Physics, 
Bharathidasan University, Tiruchirapalli - 620 024, India. 
(Dated: February 8, 2008) 

We propose a simple method to estimate the parameters involved in discrete dynamical systems 
from time series. The method is based on the concept of controlling chaos by constant feedback. 
The major advantages of the method are that it needs a minimal number of time series data and is 
applicable to dynamical systems of any dimension. The method also works extremely well even in 
the presence of noise in the time series. The method is specifically illustrated by means of logistic 
and Henon maps. 

PACS numbers: 05.45.-a,47.52.+j 
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In recent years, studies on chaotic dynamical systems 
have become extremely relevant from a physical point of 
view due to their potential applications in secure com- 
munication [1-5], cryptography [6], and so on. Also, 
much attention has been given to time series analysis 
since many physical, chemical, and biological systems 
exhibit chaotic motion in nature. The main objectives 
of time series analysis are to identify the structure of the 
equations which govern the temporal evolution of the dy- 
namical system, the number of independent variables in- 
volved, and parameters which control the dynamics of 
the system Several methods have been developed 
for modeling the dynamical systems by different authors 
ISkiSjJO, 11, l^J^J. A number of methods have also been 
proposed for estimating the system parameters based 
on the concept of synchronization [l^ lisl UtI Ii8| . 
Bayesian approach 19, 20] and least squares approach 
l21l . In this Letter, a very simple and practical method 
for estimating the system (control) parameters of discrete 
dynamical systems from the time series is developed us- 
ing the concept of controlling chaos [H HI HI . This 
method is applicable to time series obtained from a dis- 
crete system of any dimensions and can be extended to 
continuous systems also without much difficulty. The 
method can also be used for the time series which con- 
tains considerable amount of noise as well as with scalar 
time series. Further this method can be used in the field 
of controlling chaos to find the exact values of controlling 
constants (k^). 

Consider an arbitrary iV-dimensional discrete chaotic 
dynamical system (the original map), 



where k^'s are constants. The crucial idea in the con- 
struction of the modified map is that the addition of con- 
stants Ki in Eq. (1) will not affect the Jacobian of the 
original map, but it can change the original map with- 
out affecting the parameters (p) into a modified map 
exhibiting a different stable fixed point solution (other 
than the unstable fixed point of the original map). Also 
it is always possible to construct such a modified map 
by finding a suitable set of contants (^i's) which makes 
the modified map to exhibit a stable period one fixed 
point even for the parameters for which the original map 
evolves chaotically. 

Now let us start the evolution of the original and 
modified systems from a common set of initial states 
(i.e.jx'f"' = yf"''). After one time interval, the dynam- 
ics of the modified system can be represented as 
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and the dynamical variables of the original and modified 
systems can be related as 
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where c^'' = Ki. After the second interval of discrete 
time, the dynamics of the modified system can expressed 
as 

y\ ~ fi{x\ ^ + ' + C2 \ x^j^ -\- cjy^; p) -I- (5) 

and, after Taylor expansion, the relation between the 
variable of the original and modified systems becomes 
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where i = 1,2,3, ...TV, p denotes the system parameters 
of dimension M to be determined and the discrete index 
n stands for denoting the iterations. We also assume that 
the function / is sufficiently smooth. Let us construct a 
modified discrete dynamical system (the modified map) 
as 
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where 
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and X is the vector of dimension A'^. Proceeding further, 
the entire time evolution of the modified system can be 
obtained from the original system by the relation 
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and c^°^ = by our initial assumption yf'^ = xf'\ 

Next, let zf \ z"^'\ zf^~^^ be the m data set points 
of the given chaotic time series obtained for the original 
map. Then the trajectory of the modified map (which 
is constructed by adding a set of constants with the 
original map) can be obtained from the above time series 
by the relation 
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and z is a vector of dimension N. If the original system 
is one dimensional, then 



c(") =« + c("-i)^ 
dx 
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Let y* be the period one fixed point of the modified map 
obtained by Eq. (8) for the given time series data. Then 
the n*'' and (n + I)*'* iterations of the map can be ex- 
pressed as 



z'f^+c^ = y* and 



and by subtracting Eq. (10) from Eq. (11), we get 



Similarly, 
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where 0^""*"^^ c^""*"^^ and c-"^ are functions of Kj and p. 

Thus, we have obtained 2N nonlinear simultaneous alge- 
braic equations for (M + N) unknowns {N KiS and M 



TABLE I: Convergence of r and k in the logistic map 



Itera- 


Using exact time series 


Using noisy time series 


tions 


r 


K 


r 







10.00000000 


-0.50000000 


10.00000000 


-0.50000000 


1 


8.99885997 


-0.48521330 


8.99885540 


-0.48552646 


2 


7.99761227 


-0.48347134 


7.99760776 


-0.48374876 


3 


6.99636547 


-0.48126778 


6.99636104 


-0.48151029 


4 


5.99512049 


-0.47835604 


5.99511616 


-0.47856279 


5 


4.99387963 


-0.47426330 


4.99387546 


-0.47443116 


6 


3.99265086 


-0.46786433 


3.99264699 


-0.46798613 


7 


3.69848466 


-0.46349089 


3.70025299 


-0.46361426 


8 


3.67047439 


-0.46283146 


3.67212276 


-0.46295672 


9 


3.67000000 


-0.46282092 


3.67164967 


-0.46294621 


10 


3.67000000 


-0.46282092 


3.67164967 


-0.46294621 



p's), and solving them we can obtain the values of the un- 
knowns p and Ki, provided the solution exists. After esti- 
mating the unknown parameters, one can also check the 
accuracy of the estimated parameters as follows. Iterate 
the Eq. (8a) with estimated values of the parameters in 
Eq. (8b) till a fixed point solution {y*) is obtained. Then 
compare the fixed point (y*) obtained by the above iter- 
ation using the time series of the original map with the 
fixed point calculated by Eq. (2) at the estimated param- 
eters. The degree of closeness of these fixed points gives 
a measure of the accuracy in the estimated parameters. 

As an example to our method in one dimension, con- 
sider the well known logistic map 



x 



=ra;("'(l -a;(")), < a; < 1, < r < 4 (13) 

where r is the unknown system (control) parameter. 
Then the modified logistic map can be constructed as 
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where k is a constant to be determined which makes the 
modified logistic map to exhibit period one fixed point so- 
lution for the parameter where the original map exhibits 

chaotic solution. 

Let .z(2)^_ _^_j(m-i) ^jjg ^jjj^g ggj.jgg j^^g^ Q^,. 

tained from the logistic map at some arbitrary time in- 
terval for a unknown system parameter (r). Assume z^'^^ 
be the common initial state for both the original and 
modified logistic maps {i.e. x^^'' = y^^^ = z^^^) and so 
cC) = by Eq. (8a). Then after substituting three data 
points z'^^^z^'^^ and z^^^ (one can take any three succes- 
sive data) of the time series and the values of c^^\ c^^^ 
and c*^"*) calculated by making use of the Eq. (9) into the 
Eq. (12), we get 



Kr(l-K-2z(i)) 
Kr{l - 2z(2))[l -f r(l - K - 2z(i))] 
-KV[l-Fr(l-K-2z(i))]2 
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(15a) 
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The values of k and the unknown parameter (r) can be 
estimated by solving the above two nonlinear simulta- 
neous algebraic equations with an initial guess of k and 
r. For illustration purpose, we have used the numeri- 
cally generated time series of the logistic map for the 
system parameter r = 3.67 and solved the eqns. (15) by 
globlally convergent Newton's method with an ini- 
tial guess —0.5 to K and 10.0 to the parameter r. The 
convergences of the system parameter r and k are shown 



in the Table I and it shows that the estimated value r is 
3.67 which is the exact value of parameter at which the 
time series data of the logistic map is generated. 

In order to test the robustness of the method, a noisy 
time series generated by considering that the system itself 
produces some error in the data in each iteration was also 
used in the above illustration. In our analysis a noise of 
strength 10^^ is added with the eqution of the system 



TABLE II: Convergence of a, /3, ki and K2 in the Henon map 



Itera- 


Using exact time series 


Using noisy time series 


tions 


a 








a 






K2 





2.50000000 


1.50000000 


-0.10000000 


0.00000000 


2.50000000 


1.50000000 


-0.10000000 


0.00000000 


1 


2.25643981 


1.18700678 


-0.11800073 


-0.04888403 


2.25664011 


1.18682541 


-0.11787076 


-0.04876757 


2 


1.99627422 


0.88915339 


-0.14682621 


-0.10147980 


1.99662943 


0.88880258 


-0.14654089 


-0.10125436 


3 


1.72472420 


0.60769535 


-0.19746953 


-0.16839182 


1.72522662 


0.60715419 


-0.19696823 


-0.16812574 


4 


1.44784495 


0.35719041 


-0.30159310 


-0.26710065 


1.44850317 


0.35640699 


-0.30072823 


-0.26703903 


5 


1.32207233 


0.28601093 


-0.44576804 


-0.33566022 


1.31986881 


0.28336810 


-0.44645811 


-0.33730112 


6 


1.41037943 


0.30302200 


-0.45513595 


-0.33541800 


1.40939982 


0.30087213 


-0.45644206 


-0.33692901 


7 


1.40001650 


0.30000320 


-0.45914419 


-0.33391288 


1.39875604 


0.29768355 


-0.46069850 


-0.33545169 


8 


1.40000000 


0.30000000 


-0.45918942 


-0.33389223 


1.39873784 


0.29767999 


-0.46074960 


-0.33542922 


9 


1.40000000 


0.30000000 


-0.45918942 


-0.33389223 


1.39873784 


0.29767999 


-0.46074960 


-0.33542922 



that generates the time series data. For a particular set 
of data the estimated value of r is found to be 3.67164967 
after solving eqns. (15) for the same initial guess to k and 
r as before. Also, the values of parameter (r) estimated 
from the noisy time series data at various intervals (we 
have considered 1000 data points) of time is found to 
be distributed around 3.67. We have also carried out 
similar analysis for the Moran-Ricker (exponential) map 
and verified that the system parameter can be identified 
correctly both in the absence and presence of noise. 

For the illustration of our method in two dimensional 
discrete system, we consider the Henon map, 

4"+'' = l + 4"^-a(4"V, (16a) 
4"+i) = f3x["\ (16b) 

where a and /3 are control parameters to be determined, 
and the modified Henon map can be constructed as 

= l + 4")-a(y("))2 + «:i, (17a) 
4"+i) = /9yl"'+«2, (17b) 

where ki and K2 are constants which force the modified 
Henon map to exhibit period one fixed point solution for 



a set of parameters where the original Henon map shows 
chaotic behaviour. 

Let {z^^\z^°^),{z^^\z^''>),...,izt~'\zi"'-'^) be the 
data sets of time series obtained from the Henon map 
at some arbitrary interval of time for a set of unknown 
system parameters {a and /3). The starting assump- 
tion of common initial state xj""* = y^f"^ = z[^^ and 
x^°^ = y(*" = z!f^ makes cf ^ = c^"^ - by Eq. (8a). 
In the case of Henon map, the substitution of c"^\ c"^\ 
cf\ c^^ and c^'' calculated using Eq. (8b) into 

Eq. (12) leads to four nonlinear simultaneous algebraic 
equations as 

— 2q;(ki — aKi + K2 — 2aKiz[^"')z^^ 
—a{Ki — anf + K2 — 2aKiz[^'')^ 



+/3ki + K2 


= 4'^ 


^1 ' 


(18a) 


P{ki — anf + K2 — 2aKiz[^^) 


= 4' 


^2 ' 


(18b) 


2 (1) 
K2 — OiKi — 2aKiZl 


- 4'^ 


(2) 


(18c) 


I3ki 


= 4' 


(2) 


(18d) 



In this illustration, we have used the numerical time se- 
ries data of the Henon map generated for the system 
parameters a = 1.4 and [3 = 0.3 and solved the above 
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coupled (18) by the globlally convergent Newton's 
method '2^ with an initial guess a = 2.5 (3 = 1.5, 
Ki = —0.1 and K2 = 0. The convergence of a, (3, ki 
and K2 are shown in Table II and it also shows that the 
estimated values of a and (3 are 1.4 and 0.3 respectively. 
And these estimated values are in exact agreement with 
the values of the parameters for which the time series of 
the Henon map is generated. As in our previous exam- 
ple, we have solved the Eq. (18) using the time series 
data containing a random noise of strength 10~^ for the 
same initial guess. In this case, the estimated values are 
found to be a = 1.39873784 and = 0.29767999 and 
the values of parameters a and /3 estimated at various 
interval of time using the noisy data is found to be dis- 
tributed around 1.4 and 0.3, respectively. One can also 
verify that the above system parameters can be obtained 
from a scalar time series (either zi 's or Z2 's) by contruct- 
ing four algebraic equations suitably from Eq. (12) and 
making use of the system equations. For example, pa- 
rameters (a and /3) can be estimated from the scalar 
time series of zi using the four algebraic equations which 
contain zi alone in the right hand side, constructed by 
making use of c«, cf \ c['^ and cf ^ in Eq. (12). 

At this point, one may raise the question, why not 
invert directly the map (1) itself using the time series 
data so as to find the system parameters. While this 
is certainly possible in the case of exact time series, the 
extreme sensitiveness of chaotic systems to initial condi- 
tions make it an unreliable procedure in the presence of 
suitable noise. For example, in the case of logistic map 
the estimated value of r is found to be 3.7 while the orig- 
inal value is 3.67 when an 1% white noise in the range 
to 1 is introduced in the time series. On the other hand 
in our method described above, no such difficulty arises. 

Next, we wish to point out that it is possible to extend 
the analysis to identify the system itself in principle, say 
an 7V"* degree polynomial for the right hand side of Eq. 
(1). By solving sufficient number of Eqs. (12) one can 
then identify the form of the map itself. From another 
point of view, the procedure outlined here also gives a 
method to obtain the values of the controlling constants 
(Kj) for a chaotic system to a desired periodic orbit. Fi- 
nally, we have also extended the same precedure to con- 
tinuous dynamical systems for estimating the system pa- 
rameters by finding a set of differential equations which 
determine the connection between the original and mod- 
ified systems. The details will be presented elsewhere. 

To conclude, the main advantage of our method is that 
a very minimal number of time series data is sufficient 
for the accurate determination of the system parameters. 
We can check the accuracy of the estimated parameters 
by comparing the fixed point obtained by Eq. (8) using 
the time series at estimated parameters with the fixed 
point of the modified dynamical system at the same pa- 
rameters. Thus, we have developed a very simple as well 
as useful method for estimating the unknown system pa- 



rameters of the discrete dynamical systems of any dimen- 
sions and illustrated it by means of logistic and Henon 
maps. 
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